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steady-state particle current is maintained by a density gradient between 
the channel boundaries and by the influence of an external driving force. 
This leads to the development of a fluctuating but always microscopi- 
cally sharp interface between two domains of different densities which 
are fixed by the boundary chemical potentials. The internal structure of 
the interface becomes very simple for strong driving force. We calculate 
the drift velocity and diffusion coefficient of the interface in terms of the 
microscopic model parameters. 
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1 Introduction 



Diffusion of particles in long and narrow channels has a long history of the- 
oretical investigation and has recently become a focus also of experimental 
interest e.g. in the study of molecular diffusion in zeolites [I], diffusion of 
colloidal particles in confined geometry [2] or optical lattices J3j and granular 
diffusion In many such channel systems the particles cannot pass each 
other. This mutual blocking phenomenon is known as single-file effect and 
leads to subdiffusive behaviour j^HHj- Among other things the single-file ef- 
fect is responsible for low reaction effect ivity in microporous catalysts [7] and 
is thus of technical importance in chemical engineering. The single-file effect 
occurs also in biological systems, examples being the motion of ribosomes 
along the m-RNA during protein synthesis [HI E] or transport by molecular 
motors along microtubuli or actin filaments ^U]- I n single-file systems the 
longitudinal motion is the most important dynamical mode and makes such 
processes amenable to treatment by one-dimensional models [TT| I12j. 

Low- dimensional diffusive particle systems are of great interest also from a 
thermodynamic point of view. In open boundary systems, kept far from equi- 
librium by maintaining a steady state particle current, various unexpected 
kinds of critical phenomena have been discovered in recent years, includ- 
ing boundary-induced phase transitions, phase separation and spontaneous 
symmetry breaking, see ^TJ ^] and references therein for a review. These 
finite-temperature critical phenomena have no counterpart in thermal equi- 
librium since in one-dimensional systems with short range interactions there 
is no mechanism that could prevent the creation and growth of an island of 
the minority phase inside a domain of the majority phase. Therefore it is not 
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possible to have a phase-separated equilibrium state with a stable and mi- 
croscopically sharp interface between two fluctuating domains characterized 
by different values of the order parameter. 

Most of these nonequilibrium critical phenomena are not yet well-understood. 
Given the interesting diffusion properties as well as the potential for appli- 
cations to catalytic reactions it would thus be interesting to explore criti- 
cal phenomena in low-dimensional reaction-diffusion systems in more detail. 
Specifically, in this paper on one-dimensional reaction-diffusion systems we 
would like to investigate the existence and microscopic properties of inter- 
faces between coexisting nonequilibrium domains which are macroscopically 
different. 

In order to set the stage and sharpen the question we begin with some 
remarks of general nature and mention some results relevant to our approach. 
Systems of diffusing and reacting particles are usually described macroscop- 
ically by hydrodynamic equations for coarse-grained quantities like the par- 
ticle density ^3]- The density then represents the local order parameter 
specifying the spatial evolution of the macroscopic state of the system. Such 
equations are usually proposed on a phenomenological basis, paradigmatic 
examples being the Burgers equation for driven diffusive systems with par- 
ticle conservation [T3j or the Fisher equation for reactive systems without 
conservation law [To"! ITB] . These equations are in general non-linear and ex- 
hibit shocks in some cases. This means that the solution of the macroscopic 
equations may develop a discontinuity even if the initial particle density is 
smooth. This means that in these systems phase separation may occur. The 
shock represents the interface between the two thermodynamically distinct 
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phases. 

This hydro dynamic description of phase separation is, however, not fully 
satisfactory. It provides no insight into the microscopic origin of the phe- 
nomenon, and it gives no information about the internal structure of the 
shock. It could very well happen that in a particle system described on 
hydrodynamic (Eulerian) time scale by an equation which has shock solu- 
tions no corresponding microscopic discontinuity would be observable on less 
coarse-grained space or time scales which are experimentally relevant par- 
ticularly for the quasi one-dimensional systems referred to above. In order 
to understand the structure of shocks and the emergence of such nonlinear 
behaviour from the microscopic laws that govern the stochastic motion and 
interaction of particles it is therefore necessary to derive the macroscopic 
equations from the microscopic dynamics rather than postulating them on 
phenomenological grounds. 

Carrying out this programme starting from Newton's equation of motion 
constitutes a rather difficult problem. However, a substantial body of results 
of this nature has been obtained for specific one-dimensional stochastic lat- 
tice models [T7| [TH], the best-studied example being the asymmetric simple 
exclusion process (ASEP) In this basic model for a driven diffusive 

system each site k on the infinite integer lattice Z is either empty (n& = 0) or 
occupied by at most one particle (n^ = 1). A particle on site k hops randomly 
to the site k + 1 with rate D r and to the site k — 1 with rate Di, but only 
if the target site is empty. Otherwise the attempted move is rejected. The 
jumps occur independently in continuous time with an exponential waiting 
time distribution. For a single particle this is a biased random walk which on 
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large scales describes Brownian motion driven by a constant external force. 
The exclusion rule mimics a short-ranged hard-core interaction potential be- 
tween particles. In the hydrodynamic limit the system is described by the 
Burgers equation which exhibits shocks. Such a shock discontinuity may be 
viewed as the interface between stationary domains of different densities. 

Moreover, there are a number of exact results about shocks in lattice 
gas models for driven diffusive systems [2H E21 E31 EU ESI EE1 E3 EHl EH1 
EI3 EHl, in reaction-diffusion systems |2Z1 Ell E21 EH EH EH] (where shocks 
appear as Fisher waves on the macroscopic scale) and in spin-flip systems 
fI7\ EHl EH] where shocks correspond to domain walls (HHl • It has emerged 
that in all these models the macroscopic shock discontinuity originates from 
a microscopically sharp increase of the local particle density, i.e., a decrease 
of the mean distance between particles that can be observed on the scale 
of a few lattice units (which typically represents the size of particles). The 
discontinuity itself performs a biased random motion with a constant mean 
speed and diffusive mean square displacement. The existence, structure, and 
dynamical properties of microscopically sharp shocks in lattice models for 
react ion- diffusion systems are the issues on which the present work focuses. 

These results for the dynamical behaviour and microscopic properties of 
shocks have been obtained for infinite or periodic particle systems. In most 
physical applications, however, one has to study finite systems with open 
boundaries where particles are injected and extracted. This is crucial to take 
into account as - in the absence of equilibrium conditions - the boundary 
conditions determine the bulk behavior of driven systems, even to the extent 
that boundary induced phase transitions between bulk states of different 
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densities occur (371 EHl EH] • Qualitatively, the strong effect of boundary con- 
ditions on the bulk can be attributed to the presence of steady-state currents 
which carry boundary effects into the bulk of the system. Quantitatively, 
exact results for the steady state of the ASEP have helped to show that 
part of the nonequilibrium phase diagram of driven diffusive systems with 
open boundaries, viz. phase transitions of first order, can be understood 
from the diffusive motion of shocks [JOl IS], in analogy to the Zel'dovich 
theory of equilibrium kinetics of first-order transitions. In a series of recent 
papers (121 E21 HHj these considerations, originally formulated for conser- 
vative dynamics, have been extended to non- conservative reaction-diffusion 
systems with open boundaries. As in equilibrium, the nonequilibrium theory 
of boundary-induced phase transitions requires the existence of shocks which 
are microscopically sharp. Therefore, the study of the microscopic structure 
of shocks in open systems is essential for understanding boundary-induced 
first order transitions and the phase separation phenomena associated with 
it. 

After this survey we are finally in a position to precisely state the objec- 
tive of this work. All the systems studied so far allow only for the presence 
of a single species of particles. No exact results have been reported so far for 
non- stationary travelling waves in open two-component systems, i.e. where 
two diffusive particles species A, B react with each other to form an inert 
reaction product or undergo a cracking or coagulation reaction (B # 2 A). 
In order to address this question we adapt the strategy suggested in |27] to 
two-component systems: We take as initial distribution of particles a shock 
distribution with given microscopic properties and look for families of models 
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for which the shock distribution evolves into a linear combination of simi- 
lar distributions with different shock positions. Thus the information of the 
microscopic structure of the shock that one has initially is preserved for all 
times. Remarkably it will transpire that such families of reaction-diffusion 
systems exist for strong external field that drives the particles and keeps them 
in a nonequilibrium state. We remark that in a similar treatment for a dif- 
ferent family of two-component processes we have found such a phenomenon 
at some specific finite driving strength [46 . 

The paper is organized as follows: In the following section we define the 
class of models that we consider and we also define shock measures for these 
systems. In Sec. 3 we determine the families of models with travelling shock 
solutions on the finite lattice. In Sec. 4 we summarize our results and draw 
some conclusions. Some mathematical details of the calculations are given 
in the appendices. 

2 Stochastic reaction-diffusion processes 
2.1 Three-states lattice gas models 

In order to keep the physics that lead to phase-separated nonequilibrium 
states as transparent as possible we study the simplest possible setting for 
a stochastic two-component react ion- diffusion process. We consider a lattice 
gas model defined on a lattice with L sites. The state of the system at 
any given time is described by a set of "occupation numbers" n — ni, . . . , Ul 
where % = 0, 1, 2 is the local occupation number at site k. These occupation 
numbers are abstract objects and serve as mathematical labels for three 
possible local states of each lattice site. 
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The bulk stochastic dynamics are defined by nearest neighbor transitions 
between the occupation variables which occur independently and randomly 
in continuous time after an exponentially distributed waiting time. The 
mean r(n' k , n' k+1 ; n k , n k+ i) of this waiting time depends on the transition 
(nfc,rifc + i) — > (n' k ,n' k+1 ). For later convenience we introduce an integer label 

i = 3n k + n k+1 + 1 (1) 

in the range 1 < % < 9 for the occupation variables on two neighboring sites 
k and k + 1. The inverse mean transition times are the transition rates Wij. 
Here % = 3n' k + n' k+1 + 1 labels the target configuration and j is the respective 
label of the initial configuration (n k , n k+ i). We assume the bulk dynamics to 
be spatially homogeneous. The transition rates then do not explicitly depend 
on the site k. 

We require a single local conservation law where some linear function 
C (n) of the local occupation numbers is conserved under the transitions j3H] • 
It is straightforward to check that this allows for 10 nonvanishing rates Wij. 
The physical interpretation of this conservation law as charge-, mass-, or 
particle conservation respectively depends on the physical interpretation of 
the occupation numbers n k and will be given below. We present the following 
three families of models which are mathematically equivalent, but have rather 
different physical interpretations. 

2.1.1 Diffusion without exclusion 

In its most obvious interpretation the abstract occupation number n repre- 
sents the number of particles on a given site. Requiring particle conservation 
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where C(n) = n allows for 10 hopping processes with rates given as follows: 



10 # 


01 


^24, 


W 42 


20 # 


02 


W 3 7, 


W73 


12 # 


21 


™86, 


w 68 


11 # 


02 


^35, 




20 # 


11 


^57, 


W 75 . 



(2) 



Here there is no distinction between different particles, only the total number 
is recorded. Physically this process describes diffusion of a single species of 
particles in a pore system large enough to accommodate two particles in 
each pore. Thus the three states do not describe a two- component, single-file 
particle system, but a one-component system where particles can pass each 
other. This makes this process different from the previously studied two-state 
single-component systems which describe single- file diffusion |2*Hl26| l2Tj. For 
definiteness we shall focus in this paper on two-component reaction-diffusion 
systems and hence not make use of this one-component interpretation of the 
three local states. 

2.1.2 Two-species annihilation A + B # 

We define 

C{n) = 1 - n (3) 

as charge associated with the state n of a lattice. The "occupation number" 
therefore denotes an internal degree of freedom in a single-file particle sys- 
tem. The value n = corresponds to a positively charged particle (denoted 
as type A), n — 1 corresponds to a vacant site (denoted 0), and n = 2 corre- 
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sponds to occupation by a negatively charged particle (denoted as type B). 
As conservation law we require charge conservation, or, equivalently, con- 
servation of the difference of particle numbers (of positively and negatively 
charged particles). 

Since this process is mathematically equivalent to the particle conserving 
process (J2J) the dynamics of the process can be represented by the following 
ten transitions 



OA r± AO 


W 24 , 


U? 42 


BA # AB 


^37, 


W 73 


OB v± BO 


^86, 


W68 


00 # AB 


^35, 


W53 


BA # 00 


W57, 


W75 



(4) 



This is the well-studied two-component creation/annihilation process, see 
[15] for a review of some important properties and experimental significance 
of the one- dimensional pure annihilation case. The main results of this paper 
are given in terms of this process. 

2.1.3 Cracking B -> 2A 

One may switch the role of A and 0. The "occupation number" n = then 
corresponds to a vacant site 0, n = 1 corresponds to a particle of type A, 
and n = 2 corresponds to occupation by a particle of type B. We drop the 
assignment of charges to particles and instead introduce 

C(n)=n=:M (5) 
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as mass of the particles (in suitable units). A-particles thus have mass 1 and 
Z?-particles to have mass 2; the conservation law describes mass conservation. 
Under this mapping the process (JIJ) read 



AO 


^ OA 




W A2 


BO 


# OB 


W37, 


W 73 


AB ■ 


^ BA 




W68 


AA 


# OB 


^35, 




BO 


# AA 


^57, 


W75 



(6) 



The last two reactions corresponds to cracking of a molecule B with mass 
2 into two identical parts A (mass 1 each), with coagulation as reversed 
process. The third process in this list is a recombination reaction between 
neighboring reactands. 

2.1.4 Boundary conditions and continuity equation 

At the boundary sites 1,L we assume the system to be connected to some 
external reservoir with which the system can exchange particles. For defi- 
niteness we consider here and below charge conservation. The corresponding 
processes for mass conservation are obtained by changing A <-> 0. 

For injection and extraction of particles at the left boundary we introduce 
the rates : 

«2, 72, 

«3, 73, (7) 



A^O 
A^B 
O^B 
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and for the right boundary 



A^B 



fa, 







B 



S 3 , (3 3 . 



(8) 



Here and below the left rate refers to the process going from left to right, 
while the right rate is for the reversed process. The boundary rates are 
a further set of model parameters. Below we define them such that they 
are parametrized by 2 independent boundary chemical potentials which fix 
boundary densities for the conserved order parameter. 

The presence of the bulk conservation law implies a lattice continuity 
equation 



for the expectation Ck = ( C(rik) ). This quantity plays the role of a conserved 
local order parameter. The quantity jk is the current associated with the 
conservation law. It is given by the expectation of some combination of 
local occupation numbers, depending on the model under investigation, see 
below. Since we do not study here periodic systems we do not require the 
boundary sites where the system is connected to the reservoir to respect the 
conservation law. The quantities jo, 3l entering the continuity equation for 
k — 1 and k = L respectively are source terms resulting from the reservoirs. 
They are functions of the reservoir densities. The lattice continuity equation 
is the starting point for a coarse-grained hydrodynamic description of the 
time evolution of the local order parameter. 



-jjCk — jk-i — jk 



(9) 
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2.2 Master equation 

The time evolution denned above can be written in terms of a continuous- 
time master equation for the probability vector 

\P(t))=Y,P(ni,---,n L ;t)\n), (10) 

n 

where P(ni, • • • , til, t) is the distribution for the probability of finding parti- 
cles at sites 1 to L and | n ) is the basis vector in the space of configurations 
20J. The probability vector is normalized such that (s\P) = 1 with the 
summation vector ( s | = J2n ( ZL I ■ The time evolution is generated by the 
stochastic Hamiltonian H whose offdiagonal matrix elements Hn,n> are the 
negative transition rates between configurations. As required by conserva- 
tion of probability, the diagonal elements are the negative sum of transition 
rates in the respective column. 

Therefore the master equation is now described by the Schrodinger equa- 
tion in imaginary time: 

j t \P(t)} = -H\P(t)}. (11) 

with the formal solution 

\P(t))=e- m \P(0)}. (12) 

Since only nearest-neighbour interactions are included, the quantum Hamil- 
tonian H defined above has the structure 

L-l 

H = b 1 + y E h ^+i + h- (13) 

k=l 

Here b\ and bi are the boundary matrices: 

( -(«i + a 2 ) 71 
h = - a x -(7i + «3) 73 | • (14) 

13 




( 



(Si + S 2 ) 
Si 

s 2 



(31 

Wi + s 3 ) 

ft 



(h 

(P2 + (3 3 ) J L 



(15) 



The local bulk transition matrix hk,k+i with offdiagonal matrix elements — Wij 
acts non-trivially only on sites k and k + 1. Below we give hk,k+i explicitly. 

2.3 Nonequilibrium steady states 

We stress that our family of models is defined in terms of transition rates, not 
in terms of an internal energy E(n) that would determine the stationary dis- 
tribution of the process as equilibrium distribution P*(n) oc exp (—(5E(n)). 
Instead, the stationary distribution is an a priori unknown and in general 
complicated function of the transition rates. It does not in general satisfy 
detailed balance and thus represents a nonequilibrium steady state. In order 
to be able to carry out explicit computations we restrict ourselves to systems 
such that the stationary distribution of the stochastic dynamics factorizes, 
i.e., one has a product measure without correlations between the occupation 
numbers at different sites. 

Requiring the existence of a stationary product measure imposes con- 
straints both on the boundary rates and on the bulk rates. Physically, the 
constraints on the boundary rates essentially means that the chemical po- 
tentials in the two reservoirs are equal, allowing the bulk to relax into a 
current-carrying stationary state with a chemical potential determined by 
the reservoirs. In this case the origin of the current is not a gradient in 
the external chemical potential of the reservoirs, but a constant bulk driving 
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force. The conditions on the bulk rates have a less transparent and model- 
dependent physical interpretation. Once these conditions are determined 
the model is fully defined and its stationary distribution is given for equal 
chemical potentials in the reservoir. 

In the quantum Hamiltonian formalism introduced above a product mea- 
sure is given by a tensor product 

\P) = IPx)® \P 2 )®...®\P L ). (16) 

Here the three-component single-site probability vectors \Pt) has as its com- 
ponents the probabilities P{rik) of finding state n at site k. In the stationary 
distribution these probabilities are position- independent, \Pk) = \P), and the 
stationary probability vector thus has the form 

= (17) 

By definition of stationarity the stationary probability vector satisfies the 
eigenvalue equation 

H\P*) = 0. (18) 

We shall parametrize the one-site marginals P(rik) by a generalized fugacity 
z associated with the conserved quantity and an interaction parameter de- 
termined by the transition rates, see below. In formal analogy to equilibrium 
systems we shall refer to the logarithm of the fugacity as chemical potential. 

2.4 Initial conditions 

The objective of this paper is the analysis of the family of models which is 
defined by having a stationary product measure if the chemical potentials 



15 



in the reservoir are equal. However, as physical boundary conditions to be 
studied we envisage different chemical potentials in the reservoirs. The prod- 
uct measure is then no longer stationary and the questions arises what new 
properties the stationary distribution exhibits and how the system relaxes 
to its stationary distribution. Indeed, in order to avoid misunderstanding 
we stress that the product requirement on the stationary distribution with 
equal reservoir chemical potentials does not imply the absence of correlations 
during the time evolution of the more general open system with different 
reservoir chemical potentials. 

Specifically, we prepare the system initially in a state described by a 
(nonstationary) shock measure of the form 



These shock measures have single-site probabilities given by \Pi) in the left 
chain segment up to site k (chosen to match the chemical potential of the left 
reservoir) and single-site probabilities given by \P 2 ) in the remaining chain 
segment from site L — k up to site L (chosen to match the chemical potential 
at the right reservoir). 

Such a shock measure defines fully the internal structure of the shock. 
Since there are no correlations in a shock measure one may regard the lattice 
unit as the intrinsic shock width. A typical configuration has a sharp decrease 
of the mean interparticle distance across the lattice point k. In the course 
of time the measure \P(t)) = exp(—Ht)\k) changes and it is interesting 
to investigate this time evolution. For the models studied below \P(t)) 
is computed explicitly and allows for a detailed explicit calculation of all 
correlations that develop with time. 




(19) 



16 



3 Field-induced phase separation 
3.1 Stationary distribution 

Following the ideas outlined above we first search for stationary product 
solutions of the model with spatially constant single-site probabilities. By 
choosing the basis of three states as follows 




10) 



( \ 

1 

w 



\B) 




(20) 



one can conveniently write the product measure for the periodic model in 
terms of a generalized fugacity z and arbitrary constant r 



( 1 



\ rz A 



Here 



1 + Z + TZ A 



(21) 



(22) 



is the local "partition function". The quantity r parametrizes the density 
ratio of the two particle species, p B / p A = rz 2 . The fugacity z is associated 
with the conservation law, i.e., in a periodic system where the charge is 
conserved | P* ) would be stationary for any value of z. This probability 
measure is grand-canonical. The charge a = p A — p B in this ensemble has 
mean 

d n 1 — rz2 / _ _ >. 

a = \-z— In v = . (23) 

dz 



The corresponding canonical distributions with a definite value of the charge 
can be constructed in standard fashion, but we do not consider them here 
since we are dealing with an open system where the bulk fugacity is fixed 
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by the generalized chemical potentials of the reservoirs. The nonconserved 
particle density p = p A + p B in this ensemble is given by 

1 + rz 2 



P 



(24) 



The stationary distribution of the model is not known in full generality 
and we have to determine constraints on the bulk rates such that the prod- 
uct measure (}2"T]) is stationary. The transition matrix h^^+i for the bulk 
stochastic dynamics is given by 



-h 



( o 









V o 


















-(WS3 + w 73 ) 




w 73 







w 2 4 


-W24 








k,k+l — 




W35 



-(W35 + w 75 ) 


W75 





















w 37 


W57 



-(w 37 + w 57 ) 











^68 



-w 68 




(25) 



L k+1 



— n 



P*). 



(26) 



and stationarity of the product measure implies 

h k>k+1 \ P* ) = [F(n£ +1 - hi) + F'(ht 

Here F and F' are arbitrary constants and n£ and h B are number operators 
which take value 1 if there is a particle of the respective species at site k and 
zero otherwise, i.e., p A = (n A ) and p B = (n B ) independently of k due to 
homogeneity of the measure. 

In order to satisfy the relation (|27)j) for systems with open boundaries we 
can write for b\ and bi, using another arbitrary constant g 

bi \ P*} = (Fh A + F'nf + g) \ P* ) , (27) 



\ 








/ 



k,k+l 
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b L \P*) = -{Fn A L + F'nf +g)\P*). (28) 

As detailed in Appendix A one may solve the eigenvalue Eq. (|18|) and find 
F and F'\ 

F = w 2A - w i2) (29) 

F' = w m - w 68 . (30) 

Therefore the bulk rates and densities satisfy two relations due to the eigen- 
value equation (JTSJ): 

^35+^75 

r = , (31) 

W53 + w 57 

W35W57 - w 53 w 75 

W 2 4 - ^42 + ^68 - W 8 6 + W 73 - W 37 = ; . (32) 

W 35 + w 75 

The first equation (|31|) expresses the constant r in terms of the reaction 
rates. The second equation (j32j) is a constraint on the transition rates which 
we impose on the model. 

For the boundaries one needs to satisfy 

g = ^(F + F'rz 2 ). (33) 

This leaves two equations for the left boundary: 

(w A2 -w 24 )z(l + rz)-(w§ 8 -w 8 o)rz 2 + (ai+a 2 )i> -71 zv-^ 2 rz 2 v = 0, (34) 

(w 68 -w m )rz 2 (l + z)-(w i2 -w 24 )rz 2 -a 3 zv + {^ 3 + ^ 2 )rz 2 v -a 2 v = 0, (35) 
and for the right boundary one has 

(w 42 - w 24 )z(l + rz) - (w 68 - w m )rz 2 - (<5 X + b 2 )v + fay + p 2 fx = 0, (36) 

(w 68 -w m )rz 2 (l + z)-(w 42 -w 2A )rz 2 + 5 3 zv - (f3 3 + f3 2 )rz 2 v + 5 2 v = 0. (37) 
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These equations relate the boundary rates to the fugacity and moreover im- 
pose some constraints on the the boundary rates which are required for a 
proper interpretation as boundary reservoirs with fixed chemical potential. 

We remark that the given choice of nonvanishing rates is only determined 
by the conservation law and requiring stationarity of the product measure. 
Many physical processes satisfy PT-invariance, i.e., the bulk dynamics should 
be symmetric under combined time reversal and space reflection. This phys- 
ical input generalizes the equilibrium condition of detailed balance to allow 
for external driving forces. On the microscopic level of rates such driving 
forces lead to a bias in the hopping rates. Well-known examples for models 
of this kind are exclusion processes satisfying pairwise balance |Uj ■ As in the 
more general case discussed so far, the system is forced into a nonequilibrium 
steady state with a stationary current flowing in the system. Following 46 j 
we find that PT-invariance imposes the following further relations 

W75 = rw 53 , 

W35 = rw 57 . (38) 

In the calculations of the next section we do not make use of these extra 
relations. We have merely listed them for possible applications of our general 
results to specific PT-symmetric systems. 

3.2 Stationary Current and Hydrodynamics 

As remarked above the conservation law implies a lattice continuity equation 
(jDJ) for the charge current. To calculate the charge current we use the equation 
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of motion for the expected local charge density 

Jfkit) = j t [(n£) - (nf )] = j fc _! -j k . (39) 
One finds for the expected local density of A-particles 

J t ( n k) = - w ^ n l-i n k) + W42(rifc-i™°} - ^37«-i^) + w 73 (nf^n^) 
-^57«-i^fe) + WTninl_j.nl) + w 24 (n° k nf +1 ) - w i2 (n£n° k+1 ) 
+w 37 (n k 3 n k A +1 ) - w 73 {n£n^ +1 ) + w 35 (n° k n° k+1 ) - w 53 (n£n^ +1 ), 

(40) 

and for _B-particles 

j t (^k ) = ^37«_i^fc) - w^n^n^) - w m (n k ^n^) + u^faf-i"*) 
+w 35 (n° k _ 1 n° k ) - wwin^ni 3 ) - w 37 (n^n£ +1 ) + w 73 (n£n^ +1 ) 
+w 86 (n^nf +1 ) - w 68 {n^n° k+1 ) - w 57 {n^n£ +1 ) + w 75 {n° k n° k+1 ) 

(41) 

This gives the charge current 

jk = -w 24 (n° k n£ +l ) + w A2 (n£n° k+1 ) - 2w 37 (n k 3 n k A +1 ) + 2w 73 (n£n k 3 +1 ) 
-w 68 (n^ n° k+1 ) + w 86 {n° k n^ +1 ) - w 35 {n° k n° k+1 ) + Ws_{n£n% +1 ) 
-w 5 7«n^ +1 ) +w 75 (ngng +1 >. 

(42) 

In the steady state we can compute the current using the stationary distri- 
bution. One finds 

2 3 
Z 1" Z f z 

j* = (_ W24+W42 ) + [2(-w 37 +w 73 )+w 53 -w 57 }— +(w 86 -w 68 ) — +(w 75 -w 35 ) 

(43) 
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and by using (|3Tj) and the stationary condition (|32J) 

1 1 

j = ~{w i2 - w 2 i){p + a)(l - a) + -(w 8 6 - w 68 )(p - cr)(l + a), (44) 

where a and p = (n^) + (nf ) are the stationary density of charges (J2~3j) and 
particles (|2^j) respectively. 

Since the individual particle densities are not conserved the equations of 
motion for the local densities take the form 

j t (nt)=jti-j£ + S k , (45) 
j t (n%)=jti-3k+S h , (46) 

with source term 

S k = -lw 57 ({n*_X) + (nM + i)) + hrM-A) + (r« + i» 
+>35(< n *-i n 2> + " + (^nf +1 )). 

The particle currents are given by 

Jfc = - ^24«r^ +1 ) + w 42 (n^ +1 ) - w 37 (n^ nf +1 ) + w 73 (n£n% +1 ) 

1111 

- ^35(n° k n° k+1 ) + -w 53 (n^ +1 ) - -w 57 (nfn^ +1 ) + -w 75 (n° k n° k+1 ), 

(48) 



3k = w ^{ n k n k+i) - w S6 (n° k n^ +1 ) + w 37 (n^n^ +1 ) - w 73 (n^ +1 ) 

+ l w 35(n° k n° k+1 ) - -w S3 (n^ +1 ) + -w 57 (n^n^ +1 ) - -w 57 (n° k n° k+1 ). 



(47) 



(49) 
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The resulting charge current jk = j£ ~ 3k ls studied above. One may intro- 
duce also a particle current jk = jk + jjf and write 

^Pfe(^) = jk-i ~ h + 2S k . (50) 

For a coarse-grained hydrodynamic description of the time-evolution of 
the system we follow standard arguments fT7\ ITH] . We pass to a continuum 
description by making the lattice unit a (which until now had been taken 
to be a = 1) infinitesimal and we consider continuum space as x = 4. 
The coarse-grained local observables a(x, t), p(x, t) in continuous space are 
averaged over a large but finite lattice interval around the lattice point x 
and therefore given by the expected local densities a x (t), p x (t). We consider 
Eulerian scaling t' = ta with rescaled macroscopic time if. In the continuum 
limit the two equations for a and p then take the form (to leading order in 
the lattice constant a) 

d v a(x,t') = -d x j(a,p), (51) 

d t ,p(x, t') = -d x j(a, p) + R(a, p)/a + R(a, p), (52) 
where because of local stationarity 

R(cr,p) = -~(w 5 7 + w 53 )(p + a)(p - a) + 2(w 75 + w 35 )(l - p) 2 , (53) 

R{a,p) = {w 57 - w 53 )^[{p + a)d x {p - a) - (p - a)d x (p + a)}. (54) 

The space-time dependence of R and R is implicit in arguments a(x, t'),p(x, if) 
In this limit, when time and space are large, the term contained R(a, p) 
in the equation for p becomes large enough to make the two other terms neg- 
ligible. Therefore p(x,t) reaches its stationary state very fast, in agreement 
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with the notion that non-conserved local degrees of freedom have attained 
their stationary values under hydrodynamic scaling. This implies that in the 
stationary state R(cr, p) = 0, from which we obtain the stationary particle 
density 



for a given value of charge a. Therefore p takes at any instant of (macro- 
scopic) time a special value p* which is a function of a. The remaining slow 
dynamical mode is the charge, the evolution of which is thus governed by the 
hydrodynamic equation 



In the second equation j(p*, u) is the stationary current (jH| . This evolution 
equation is a nonlinear partial differential equation which can be solved by 
the method of characteristics. Because of the nonlinearity the solution may 
develop shocks in the charge distribution and we now turn to the investigation 
of these shocks on microscopic scale. 

3.3 Shock measures 

We assume that the initial distribution of charges exhibits a shock which 
on microscopic scale is represented by a shock measure (see Figure Q). We 
represent a shock measure with a shock in the fugacities between sites k and 
k + 1 as 



In this model with open boundary condition, the first (second) fugacity 
matches the fugacity in left (right) boundary. 




(55) 




(56) 




(57) 
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Figure 1: Coarse grained density profiles of a shock measure with shock 
between sites k,k + 1. 



Now we investigate the possibility that in analogy to the processes con- 
sidered in |2*7| |35J |2H] the family of shock measures | k ) is closed under the 
time evolution t. This means that the initial measure evolves into a lin- 
ear combination of shock measures after time t. This condition requires H 
which generates the time evolution to satisfy the following equation after an 
infinitesimal step 

^\k) = d 1 \k-l) + d 2 \k+l)-(d l + d 2 )\k). (58) 

We remark that this equation for the full particle distribution is mathemat- 
ically equivalent to the evolution equation of a single-particle random walk 
with hopping rate di to the left and d 2 to the right. Thus, if (J5HJ) can be sat- 
isfied, the shock in the initial distribution remains microscopically preserved 
at all times, but its position performs a random walk with shock hopping 
rates di to the left and d 2 to the right respectively. 
For further analysis we define 

h i>i+ i = h i>i+1 + F(nf - nf +1 ) + F\nf - nf +1 ), (59) 

&! = &!- Fnf - F'nf , (60) 
b L = b L +Fn£ + F'nl (61) 
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Using 

hi,*i\k) = for i^k, (62) 
=9i\ k), b L \k) = -g 2 \k). (63) 

with 

9l = -F--F' r -^, (64) 

g 2 = -F--F' r A (65) 
v 2 v-i 

yields 

-H\k) = -(J2hi+i + ~bi + bL)\k) = (-h ktk+1 -gi + g 2 )\k). (66) 

i 

Together with (jSHj) we thus find 

{-h k ,k+\ + d 1 + d 2 - gi + g 2 )\k) - d±\ k - 1) - d 2 \ k + 1 ) = 0. (67) 

The quantities gi^ are obtained from the boundary conditions (Appendix 
B). 

This is a set 9 equations for the bulk rates. We have found a solution (see 
Appendix B) with w 2 4 = 0. Putting this into the 9 equations (|B-1J) - (|B-9J) 
one finds after some straightforward algebra 

d 2 = z 2 = 0, (68) 

d 1 = * = ™« (69) 
z/i z/i 

^57 = w Z7 = 0, (70) 

W 8 6 = W 68 . (71) 

In this model there is a strong driving force for the positive particles that 
leads them to move only to the right as in the totally asymmetric simple 
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exclusion process, z 2 = means that in the right branch of the shock the 
lattice is completely filled with positive particles (see Fig. |2J). Hence incoming 
A-particles which react with B particles in the left branch of the shock hit 
the pure A-domain where they stop because of the single-file (exclusion) 
condition. The shock that separates the two domain moves only to the left 
with rate d\. Hence its mean velocity v s and diffusion coefficient D s are 
determined by the density and hopping rate only of the A-particles in the 
left domain 

v s = 2D S = w 42 pf. (72) 





Pa 


o {1) 
Pb 




Pa 


(2) 

Pb 





1 k L 

Figure 2: Density profile of a shock measure in the case z 2 = 



The interpretation of this result for the cracking process is readily avail- 
able by interchanging the role of positive particles and vacancies. The right 
branch of the shock is the empty lattice where no reactions are going on. The 
left branch is active. All particles are driven to the left so that the inactive 
region grows diffusively with drift and fluctuations determined by (|T2*|) and 
Pi replaced by the vacancy density in the active domain. 

We note that PT-invariance of the special model with u>24 = leads to 

^35 = 0, (73) 
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and to the stationary state condition 

W42 = W 73 +W 53 . (74) 

The properties of the shock are not effected by PT-invariance. 

4 Conclusions 

We have studied the dynamics of a family of one- dimensional driven two- 
component reaction-diffusion processes with open boundaries on microscopic 
lattice scale and derived a hydrodynamic description on coarse grained Eu- 
lerian scale. This is the first main result, see Eqs. (JUJ), ([55]). (foT)j) . The 
hydrodynamic equation is nonlinear and therefore admits shock solutions, 
corresponding to phase-separated states of the system. This generalizes one- 
dimensional field-induced phase separation that has been studied in some 
detail for lattice fluids in thermal equilibrium |49j . 

As the second main result we have obtained for a subset of models with 
very strong driving force detailed knowledge about the microscopic structure 
of the shock. The transition between the two phases is maximally sharp on 
lattice scale and the shock position performs a biased random walk with drift 
velocity and diffusion coefficient ()72|). Therefore, as observed in other models, 
shocks behave like collective single-particle excitations already on the lattice 
scale - not only after coarse-graining where all the microscopic features of the 
shock are lost. To our knowledge together with these are the first results 
of this nature obtained for two-component reaction-diffusion processes. The 
mapping to different models, in particular to the partial exclusion process, 
suggests that this feature is not specific to single-file diffusion. The driving 
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force that is required to produce a maximally sharp interface depends on the 
model under consideration 

For a more general choice of rates the microscopic time evolution of the 
shock structure is more complicated. By analogy with general considerations 
of microscopic shock stability j|20j and exact results for the asymmetric simple 
exclusion process [22] it is natural to expect a microscopically sharp interface, 
but with some extended structure that depends on the strength of the driving 
force. It would be very interesting to investigate experimentally in effectively 
one-dimensional driven reaction-diffusion systems with a conservation law 
the existence of field-induced phase separation and force-dependence of the 
internal width of the domain boundary. 
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Appendix A: Stationarity condition 

Assuming product measure as stationary solution, we have 





With (|59 |) -([61 |) . where nf and nf are number operators 



/ 1 \ / 

nf = 000 , n f = 00 

V o o a). \o o 





and eigenvalue equation 



H\P*) = 0, 



(A-3) 
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we write 

h i ,i+i\P*) = (b 1 + b L )\P*) = 0. (A-4) 
h iti+ i in terms of arbitrary constants F and F' is given by 



( 
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(A-5) 

where 

9x = ^-^' + ^53 + ^73, (A-6) 

2 = -F + F' + w 37 + w n7 . (A-7) 

Substituting /i^+i in the Eq. (A-7) yields 5 independent equations. One gets 
F and F' by solving following equations 

(F + w 42 - w 24 )z = 0, (A-8) 
(-F' + w S e - w 68 )rz 2 = 0. (A-9) 



Hence F and F' are 



F = w 24 - w 42 , (A-10) 
F' = w 86 -w 68 . (A-ll) 
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Three remained equations which have to be satisfied are 

(w 37 - Qi)rz 2 + w 35 z 2 = 0, (A-12) 
{w 73 - 6 2 )rz 2 + w 75 z 2 = 0, (A-13) 
O53 + w 57 )rz 2 - (w 35 + w 75 )z 2 = 0. (A-14) 

From Eq. (A- 17) we obtain 

r = W35 + W75 . (A-15) 
W53 + w 57 

Subtracting Eq. (A-15) from Eq. (A-16) yields second stationary state con- 
dition 

- ^53^75 /A 1 r>\ 

W24 - W42 + w 6S - ^86 + w 73 - w 37 = , (A-16) 

W35 + w 75 

where the sum of (A-15) and (A-16) is already satisfied. 

This model is assumed to have open boundaries, therefore b\ and bi in 
terms of injection and extraction rates are given by Eq. (JT3J) and Eq. ()15j) . 
For satisfying Eq. (A-6) for the model with open boundaries one writes 

h\P*) = (Fnf + F'nf + g)\P*), (A-17) 

b L \ P* > = -{Fh A L + F'h B L +g)\P* >. (A-18) 

where g is an arbitrary constant. Eq. (A-22) for the left boundary leads to 
three equations 

(«i + a 2 ) - jiZ - ^rz 2 = F + g, (A-19) 
-ai + (71 + a 3 )z - 7 3 r2; 2 = gz, (A-20) 
-a 2 - a 3 z + ( 72 + l3 )rz 2 = (F' + g)rz 2 . (A-21) 
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We then obtain g 



F + F'rz 2 



9 = 



v 



act 

71 + "3 ^73- 



(A-22) 



One also can obtain two conditions for boundary rates, which for the left one 

(w4 2 -W24)z(^+rz)-(w 6 8-Ws6)rz 2 + (ai + a 2 )is-'yizi'-'y2rz 2 i' = 0, (A-23) 

(wQs-wsQ)rz 2 (l+z)-(w i 2-W2A)rz 2 -a3Zv+(^3+^2yz 2 v-ot2i' = 0, (A-24) 
and for the right boundary 

(w A2 - w 2 a)z{1 + rz) - (u>68 - ws G )rz 2 - (5 1 + 8 2 )v + j3\v + (3 2 fi = 0, (A-25) 

(w es - w 8e )fi(l - n) - O42 - 1024)^ - {Jh + + 73^ + 5 2 v = 0. (A-26) 

Appendix B: Random walk conditions for the 
shock 

Explicitly the equations (J57)l that solve the random-walk condition for the 
shock are given by 



(S - w 35 - w 75 )z 1 z 2 + w 53 rz% + w 57 rz 2 - d x —z\ - d 2 —z 2 = 0, (B-2) 




(B-l) 



v 2 v x 




(B-3) 




(B-4) 



(B-5) 
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^2 jy-^ 

(S - w 68 )z 2 + w &8 z x - di — d 2 — — = 0, (B-6) 

Z\ v 2 z 2 v 2 

(S - wgejzx + w 86 z 2 - d\— d 2 — — = 0, (B-7) 

Z\ v 2 z 2 v 2 

(S - w 37 - w 53 - A)rz\ + w 35 z l z 2 + w 37 rz\ - d x rz\ — - d 2 rzf — = 0, (B-8) 

V 2 V\ 

(S - w 73 - w 57 + A)rz\ + w 75 ziz 2 + w 73 rz\ - dirz 2 — - d 2 rz 2 — = 0, (B-9) 

V 2 V\ 

where for compact notation we have introduced 

S = di + d 2 - g x + g 2 ; A = w 24 - w A2 + w 68 - w 86 + w 73 - w 37 . (B-10) 

These relations can be rewritten as 4 independent relations between the 
hopping rates and the fugacities 

w 2i = w 68 =p, (B-ll) 

^42 = w m = q, (B-12) 

V - = 4 = X 2 , (B-13) 

q z{ 

S = p + q. (B-14) 
and two equations for the shock hopping rates 

dt = q—, (B-15) 
d 2 = p— . (B-16) 

^2 

To be more specific, solving Eq. (jB-l|) and ()B-4j) - (jB-5j) yields Eq. ()B-15|) and 
Eq. (lB^T6l) for d x and d 2 , from these two and Eq. (|B^3|) and (|B^6|) - (|B^7j) . we 
obtain ()B-13|) . a relation between rates and densities, then using Eq. (jB-l|) 
with above results yields Eq. (|B7TT|l . 1B^T2|) and (|B^T4l . 
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Using these 6 relations f|B^TT|l-(fB^T6|l. equations (jE^jl . (|rT9|) and Eq. (lH^2l) 

respectively lead to the following relations for the so far undetermined rates 

(p - q)(l - + u, 35 (^| - 1) = 0, (B-17) 

(? - P)(l " + ^ 75 (^ - 1) = 0, (B-18) 

(p + g ) r -i + tUfi3 ( _ 1) + W57 ( jl _ 1) = 0. (B-19) 
V 9 VP 

Simplifying Eq. (|B-14j) by using Eqs. (J53)l - (J55Jl for #1 an d Qi yields following 
more explicit relation between r and X 

r = (lW' (B ' 20) 

This relation on r together with Eqs. ()B-17|) - (jB-19|) and the stationary state 
equation (|32p. implies that X = 0. This solved by p = z 2 = 0. 
The boundary equations (|64| ) -()65| ) lead to 

P-Q f 2 n 
= -«i 73^1 + (71 + a 3 ), (B-21) 



92 = {rz 2 - 1) 

^2 

= 5i- + ^ 2 -(/?i + 5 2 ), (B-22) 

for ^! and #2 respectively. 
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